(************************************************************************)
(*                                                                      *)
(*    Release 3.00                                                      *)
(*                                                                      *)
(*    CoreMATH - Common math support module.                            *)
(*                                                                      *)
(*    Copyright (C) 1989, 1990 Jensen and Partners. All Rights Reserved *)
(*    Written by Simon Knight.                                          *)
(*                                                                      *)
(*    Copyright (C) 1991 Psion PLC. All Rights Reserved                 *)
(*    Modified by NSM.                                                  *)
(*                                                                      *)
(************************************************************************)

include "corelib.inc"
include "epocdefs.inc"

module CoreMath

(****************************************************************************)
segment _BSS(BSS,28H)

public _matherr :   org CodePtrSize
public __SH_sword : org 2
public __exc_addr : org 2

segment _BSS(BSS,28H)

public __fac :      org 10

segment _DATA(DATA,28H)

public _HUGE : db 098H,060H,0B9H,0D7H,0FFH,0FFH,0EFH,07FH

segment _DATA(DATA,28H)

public _LHUGE : db 0FFH,0FFH,0FFH,0FFH,0FFH,0FFH,0FFH,0FFH,0FEH,07FH

(****************************************************************************)
section

extrn FIDRQQ
extrn _matherr
extrn FIWRQQ
extrn __seterrno


segment _TEXT(CODE,00048H) (*segment length = 000A2H*)

public __report_math_error :
    push  bp
    mov   bp,sp
    sub   sp,30H
    push  di
    push  si
    push  cx
    fwait
    fld     st(0),st(1)
    fwait
    fstp    tbyte [bp][-30H],st(0)
    fwait
    fst     st(1),st(0)
    mov   cx,ax
    cmp   word [_matherr],0
    je	  $1
    xor   si,si
$0:
    mov   di,si
    add   di,bx
    mov   al,cs:[di]
    mov   [bp][si][-0AH],al
    inc   si
    cmp   si,0AH
    jl	  $0
    mov   [bp][-26H],cx
    fwait
    fld     st(0),st(5)
    fwait
    fstp    qword [bp][-22H],st(0)
    nop
    fwait
    fwait
    fld     st(0),st(6)
    fwait
    fstp    qword [bp][-1AH],st(0)
    nop
    fwait
    lea   ax,[bp][-0AH]
    mov   [bp][-24H],ax
    lea   ax,[bp][-26H]
    call  word [_matherr]
    test  ax,ax
    jne   $3
$1:
    cmp   cx,1
    jne   $4
    mov   ax,EDOM
$2:
    call  __seterrno
    fwait
    fld     st(0),st(1)
    fwait
    fstp    st(1),st(0)
    jmp   $5
$3:
    fwait
    fld     qword st(0),[bp][-12H]
    fwait
    fstp    st(1),st(0)
    jmp   $5
$4:
    mov   ax,ERANGE
    jmp   $2
$5:
    fwait
    fld     tbyte st(0),[bp][-30H]
    fwait
    fstp    st(2),st(0)
    pop   cx
    pop   si
    pop   di
    mov   sp,bp
    pop   bp
    ret 0

segment _TEXT(CODE,28H)

two: dw 2

extrn __report_math_error

__lib_error:
  push bp
  push ax
  push bx
  push cx
  push dx
  mov bp,sp
  sub sp,2
  test byte cs:[si][4], f_st5 (* is arg1 in st5? *)
  jz not_st5
  fld st(0),st(5)
  fstp st(1),st(0)
  jmp st5_done
not_st5:
  fst st(5), st(0)
st5_done:
  ftst st(0)
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 1
  mov ax, cs:[si][2]
  jz me_positive
  mov ax, cs:[si][0]
me_positive:
  call ax
  fstp st(1),st(0)
  or ax,ax
  jz me_ok
  test byte cs:[si][4], f_two
  jnz has_two
  fldz st(0)
  fstp st(7),st(0)
has_two:
  lea bx,[si][5]
  call __report_math_error
me_ok:
  mov sp,bp
  pop dx
  pop cx
  pop bx
  pop ax
  pop bp
  ret 0

segment _TEXT(CODE,28H)
err_domain:
  fldz st(0)
  mov  ax, _DOMAIN
  ret 0

segment _TEXT(CODE,28H)
err_underflow:
  fldz st(0)
  mov ax, _UNDERFLOW
  ret 0

segment _TEXT(CODE,28H)
err_plus_infinity:
  extrn _HUGE
  fld qword st(0), [_HUGE]
  mov ax, _OVERFLOW
  ret 0

segment _TEXT(CODE,28H)
err_minus_infinity:
  call err_plus_infinity
  fchs st(0)
  ret 0

segment _TEXT(CODE,28H)
one:dw 1
check_abs_one:
  fld st(0),st(0)
  fabs st(0)
  ficomp word st(0), cs:[one]
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 40H
  jz ab1_err
  ret 0
ab1_err:
  add sp,2
  jmp err_domain

segment _TEXT(CODE,28H)
ldexp_err_pos:
  fld st(0),st(6)
  ftst st(0)
  fstsw [bp][-2]
  fstp st(0), st(0)
  test byte [bp][-1], 1
  jnz err_underflow
  jmp err_plus_infinity
ldexp_err_neg:
  fld st(0),st(6)
  ftst st(0)
  fstsw [bp][-2]
  fstp st(0), st(0)
  test byte [bp][-1], 1
  jnz err_underflow
  jmp err_minus_infinity

segment _TEXT(CODE,28H)
public _hypotl:
public _hypot:
hypotCommon:
  fmul st(0), st(0)
  fld  st(0), st(6)
  fmul st(0), st(0)
  faddp st(1), st(0)
  fsqrt st(0)
  ret	0

segment _TEXT(CODE,28H)

public _ceill:
public _ceil:
  push bp
  mov bp,sp
  sub sp,2

ceilCommon:
  fstcw [bp][-2]
  fwait
  and byte [bp][-1], ~4   (* round towards +inf *)
  fldcw [bp][-2]
  frndint st(0)
  or byte [bp][-1], 4     (* restore to chop *)
  fldcw [bp][-2]
  fwait
  mov sp,bp
  pop bp
  ret	0

segment _TEXT(CODE,28H)
public _floorl:
public _floor:
  push bp
  mov bp,sp
  sub sp,2
floorCommon:
  fstcw [bp][-2]
  fwait
  and byte [bp][-1], ~8  (* round towards -inf *)
  fldcw [bp][-2]
  frndint st(0)
  or byte [bp][-1], 8    (* restore to chop *)
  fldcw [bp][-2]
  fwait
  mov sp,bp
  pop bp
  ret	0

segment _TEXT(CODE,28H)
public __round:
  frndint st(0)
  ret 0

segment _TEXT(CODE,28H)
public __bcd:
  push bp        (* new *)
  mov bp, sp     (* new *)
  sub sp, 2      (* new *)
  fstcw [bp][-2] (* new *)
  (* ensure round to nearest (80387 compatibility problem) *)
  fwait          (* new *)
  push [bp][-2]  (* new *)
  and word [bp][-2], ~0C00H  (* new *)
  fldcw [bp][-2]             (* new *)

  fld st(0), st(0)
  xchg ax,bp
  fbstp  tbyte [bp], st(0)
  xchg ax,bp

  pop [bp][-2]   (* new *)
  fldcw [bp][-2] (* new *)
  fwait          (* new *)
  mov sp, bp     (* new *)
  pop bp         (* new *)
  ret 0

segment _TEXT(CODE,28H)
(* compute 2^st(0) , result is represented as st(0) * 2^st(4) *)
four: dw 4
pow2:
  fld1 st(0)
  fxch st(0), st(1)
  fst  st(5), st(0)
  frndint st(0)
  fsub  st(0), st(1)
  fxch  st(0), st(5)
  fsub  st(0), st(5)  (* we now know that st(0) is in range [0..2) *)
  fidiv word st(0), cs:[four] (* now in range [0..0.5] *)
  f2xm1 st(0)
  faddp st(1), st(0)
  fmul  st(0), st(0)
  fmul  st(0), st(0)
  ret 0

segment _TEXT(CODE,28H)
public _fmodl:
public _fmod:
  push bp
  mov bp,sp
fmodCommon:
  sub sp,2
  fld st(0), st(6)
  ftst st(0)
  fstsw [bp][-2]
  test byte [bp][-1],40H
  jnz fmod_zero
  fxch st(0), st(1)
fmod_more:
  fprem st(0)
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 4
  jnz fmod_more
fmod_zero:
  fstp st(1),st(0)
  mov sp,bp
  pop bp
  ret	0

segment _TEXT(CODE,28H)
public _frexpl :
public _frexp :
  push bp
  mov bp,sp
  sub sp,2
  xchg bx,ax
frexpCommon:
  ftst st(0)
  fstsw [bp][-2]
  fwait
  test byte [bp][-1],40H
  jnz frexp_zero
  fxtract st(0)
  fxch   st(0), st(1)
  fistp  word [bx], st(0)
  fwait
  inc    word [bx]
  fidiv word st(0), cs:[two]
frexp_exit:
  xchg bx,ax
  mov sp,bp
  pop bp
  ret	0
frexp_zero:
  and    word [bx], 0
  jmp frexp_exit

segment _TEXT(CODE,28H)
public __clear87 :
  push  bx
  push  ax
  mov   bx, sp               (* get current status word into ax *)
  fstsw ss:[bx]
  fclex
  pop   ax
  pop   bx                   (* clear status word *)
  ret	0

segment _TEXT(CODE,28H)
public __fpreset :
  ret	0

segment _TEXT(CODE,28H)
public __status87 :
  push bx
  push  ax
  mov   bx, sp               (* get current status word into ax *)
  fstsw ss:[bx]
  fwait
  pop   ax
  pop   bx
  ret	0

segment _TEXT(CODE,28H)
public  __control87 :
new  = -2
mask = -4
locals = 4

  push  bp
  mov   bp, sp               (* get current control word into ax *)
  sub   sp, locals
  push bx
  push cx
  push dx
  mov  [bp][new], ax
  mov  [bp][mask], bx
  push  ax
  mov   bx, sp
  fstcw ss:[bx]
  fwait
  pop   ax                   (* ax has original value *)

  mov   cx, [bp][mask]
  jcxz  con87Done            (* mask is zero - return *)
  mov   dx, [bp][new]
  and   dx, cx               (* ctrl = ((ctrl &~mask) | (new & mask)) *)
  not   cx
  and   ax, cx
  or    ax, dx
  mov   [bp][new], ax
  fldcw [bp][new]
con87Done:
  pop  dx
  pop  cx
  pop  bx
  mov  sp, bp
  pop bp
  ret	0

segment _TEXT(CODE,28H)
Pi_By_4: dw 0C235H, 02168H, 0DAA2H, 0C90FH, 03FFEH
public _sin :
public _sinl :
  push  bp
  mov   bp, sp
  sub   sp, 2
  push  cx
  mov  cl, 0
  jmp  trig
public _cos:
public _cosl :
  push  bp
  mov   bp, sp
  sub   sp, 2
  push  cx
  mov  cl, 2
  jmp  trig
public _tan :
public _tanl :
  push  bp
  mov   bp, sp
  sub   sp, 2
  push  cx
  mov  cl, 4
trig:
  push  ax
trigCommon:
  fxam  st(0)    (* '87 works while we enter routine*)
  fstsw  [bp][-2]
  fwait
  mov  ah, [bp][-1]
  fabs st(0)
  fld tbyte st(0), cs:[Pi_By_4]
  fxch st(0),st(1)
  fprem st(0)
  mov  ch, ah
  and  ch, 2
  shr  ch, 1  (* ch := sign *)
  fstsw  [bp][-2]
  fwait
  mov  ax, [bp][-2]
  mov  al, ah
  and  al, 3    (* select the C3, C1, and C0 status bits*)
  shl  ah, 1
  shl  ah, 1
  rcl  al, 1    (* arrange them in the ls 3 bits of BX*)
  add  al, 0FCH
  rcl  al, 1    (* al now contains the octant number*)
  cmp  cl, 2      (* are we doing Cosine ?*)
  jne  trig_notCosine
  add  al, cl     (* Cos (x) = Sin (x + 2 octants)*)
  mov  ch, 0      (* Cosine has even symmetry around 0*)
trig_notCosine:
  and  al, 7
  test  al, 1
  jz  trig_evens
  fsubp st(1),st(0)  (* overwrites pi/4 in st(1), pops.*)
  jmp    trig_ptan
trig_evens:
  fstp  st(1),st(0)
trig_ptan:
  fptan st(0)
  cmp  cl, 4
  je  trig_tangent
  test  al, 3
  jpe  trig_sineOctant
trig_cosOctant:
  fxch st(0), st(1)
trig_sineOctant:
  fmul  st(0), st(0)
  fstp  st(7), st(0)
  fld   st(0), st(0)
  fmul  st(0), st(0)
  fadd  st(0), st(7)
  fsqrt st(0)
  shr  al, 1
  shr  al, 1
  xor  al, ch
  jz  trig_sineDiv
  fchs st(0)
trig_sineDiv:
  fdivp st(1), st(0)
  jmp trig_end
trig_tangent:
  mov  ah, al
  shr  ah, 1
  and  ah, 1
  xor  ah, ch
  jz  trig_tanSigned
  fchs st(0)
trig_tanSigned:
  test  al, 3
  jpe  trig_ratio
  fdivrp st(1),st(0)
  jmp trig_end
trig_ratio:
  fdivp st(1), st(0)
trig_end:
  pop ax
  pop cx
  mov sp, bp
  pop bp
  ret	0

segment _TEXT(CODE,28H)
Pi_By_2: dw 0C235H, 02168H, 0DAA2H, 0C90FH, 03FFFH

segment _TEXT(CODE,28H)
public _atan2l:
public _atan2:
  push bp
  mov bp,sp
  sub sp,4
atan2Common:
  call atan20
  mov sp,bp
  pop bp
  ret	0
atan20:
  ftst  st(0)
  fstsw [bp][-2]
  fwait
  test  byte [bp][-1],1
  je    atan21
  fchs  st(0)
  call  atan21
  fchs    st(0)
  ret 0
atan21:
  fld   st(0), st(6)
  ftst st(0)
  fstsw [bp][-4]
  fwait
  test  byte [bp][-3],1
  je    atan22
  fchs  st(0)
  call  atan22
  fldpi st(0)
  fsubrp st(1),st(0)
  ret 0
atan22:
  test byte [bp][-3],40H
  jnz maybe_zero_zero
not_zero_zero:
  fcom  st(0), st(1)
  fstsw [bp][-2]
  fwait
  test  byte [bp][-1],1
  jz    atan23
  fxch   st(0), st(1)
  fpatan st(1), st(0)
  fld    tbyte st(0),cs:[Pi_By_2]
  fsubrp st(1),st(0)
  ret 0
atan23:
  fpatan st(1),st(0)
  ret 0
maybe_zero_zero:
  test byte [bp][-1],40H
  jz not_zero_zero
  fstp st(0), st(0)
  push si
  mov si, atan2_err_data
  call __lib_error
  pop si
  ret 0
atan2_err_data: dw err_domain, err_domain; db 0, "atan2", 0

segment _TEXT(CODE,28H)

public _atanl:
public _atan:
  push bp
  mov bp,sp
  sub sp,2
atanCommon:
  call atan0
  mov sp,bp
  pop bp
  ret 0
atan0:
  ftst  st(0)
  fstsw [bp][-2]
  fwait
  test  byte [bp][-1],1
  je    atan1
  fchs  st(0)
  call  atan1
  fchs    st(0)
  ret 0
atan1:
  fld1  st(0)
  fcom  st(0), st(1)
  fstsw [bp][-2]
  fwait
  test  byte [bp][-1],1
  jnz    atan2
  fpatan st(1), st(0)
  ret 0
atan2:
  fxch   st(0), st(1)
  fpatan st(1),st(0)
  fld    tbyte st(0),cs:[Pi_By_2]
  fsubrp st(1),st(0)
  ret 0

segment _TEXT(CODE,28H)
public _logl :
public _log :
  push bp
  mov bp,sp
  sub sp,2
logCommon:
  ftst st(0)
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 41H
  jnz log_error
  fldln2 st(0)
  fxch st(0),st(1)
  fyl2x st(1),st(0)
log_exit:
  mov sp,bp
  pop bp
  ret 0
log_error:
  push si
  mov si, log_err_data
  call __lib_error
  pop si
  jmp log_exit
log_err_data: dw err_domain, err_minus_infinity; db 0, "log", 0

segment _TEXT(CODE,28H)
public _log10l :
public _log10 :
  push bp
  mov bp,sp
  sub sp,2
log10Common:
  ftst st(0)
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 41H
  jnz log10_error
  fldlg2 st(0)
  fxch st(0),st(1)
  fyl2x st(1),st(0)
log10_exit:
  mov sp,bp
  pop bp
  ret 0
log10_error:
  push si
  mov si, log10_err_data
  call __lib_error
  pop si
  jmp log10_exit
log10_err_data: dw err_domain, err_minus_infinity; db 0, "log10", 0

segment _TEXT(CODE,28H)
__math_magic_start:

(* routines which use 8087 exception handling *)
(* to avoid explicit checks *)

segment _TEXT(CODE,28H)
public _fabsl :
public _fabs:
  mov  byte ss:[InMathLib], 1
fabsCommon:
  fabs st(0)
  mov  byte ss:[InMathLib], 0DEH
  ret  0

segment _TEXT(CODE,28H)
public _sqrtl :
public _sqrt :
  mov  byte ss:[InMathLib], 1
sqrtCommon:
(* Return the square root of x, which must be greater than zero. *)
  fsqrt st(0)   (*invalid op possible*)
  fwait
  mov  byte ss:[InMathLib], 0DEH
sqrt_exit:
  ret 0
  dw magic, sqrt_exit, err_domain, err_domain; db 0, "sqrt", 0

segment _TEXT(CODE,28H)
public _expl :
public _exp :
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,2
expCommon:
  fst st(5), st(0)            (* save parameter for error processing *)
  fldl2e st(0)
  fmulp st(1), st(0)
  fist word [bp][-2], st(0) (* to check that it is not too large *)
  fwait
  call pow2
  fld st(0), st(4) (* retrieve integral part *)
  fxch  st(1), st(0)
  fscale st(0), st(1)
  fstp st(1),st(0)
exp_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, exp_exit, err_underflow, err_plus_infinity; db f_st5, "exp", 0

segment _TEXT(CODE,28H)
public __pow_int:
  push  cx
  push  bx
  mov   bx,ax
  call  powi0
  pop   bx
  pop   cx
  ret   0

powi0:
  test  ax,ax
  jg    powi3
  jne   powi2
  fld1  st(0)
  fstp  st(1),st(0)
  ret 0

powi2:
  neg   ax
  call  powi3
  fld1    st(0)
  fdivrp  st(1),st(0)
  ret 0

powi3:
  mov   cx,17
  fld   st(0), st(0)
powi4:
  dec   cx
  add   ax,ax
  jnc   powi4
  jmp   powi7
powi6:
  fmul  st(0),st(0)
  add   ax,ax
  jnc   powi7
  fmul  st(0),st(1)
powi7:
  loop  powi6
  fstp  st(1),st(0)
powi7_exit:
  ret 0
  dw magic, powi7_exit, powi_err_neg, powi_err_pos; db f_st5, "pow", 0
powi_err_neg:
  test bx,bx
  jl err_underflow
  test bx,1
  jz err_plus_infinity
  jmp err_minus_infinity
powi_err_pos:
  test bx,bx
  jl err_underflow
  jmp err_plus_infinity

segment _TEXT(CODE,28H)
public _powl :
public _pow :
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,4
powCommon:
  ftst st(0)
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 41H
  jnz pow_integral
  fst st(5), st(0)            (* save parameter for error processing *)
  fld1 st(0)
  fxch st(0),st(1)
  fyl2x st(1),st(0)
  fmul st(0), st(6)
  fist word [bp][-2], st(0) (* to check that it is not too large *)
  fwait
pow1:
  call pow2
  fld   st(0), st(4) (* retrieve integral part *)
  fxch  st(1), st(0)
  fscale st(0), st(1)
  fstp st(1),st(0)
pow_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, pow_exit, err_pow_neg, err_pow_pos; db f_st5+f_two, "pow", 0

pow_integral:
  test byte [bp][-1], 40H
  jnz pow_dom_err
  fld st(0), st(6)
  fist word [bp][-4], st(0)
  ficomp  word st(0), [bp][-4]
  fstsw [bp][-2]
  fwait
  test byte [bp][-1], 40H
  jz pow_dom_err
  mov ax, [bp][-4]
  call __pow_int
  jmp pow_exit
  dw magic, pow_exit, err_domain, err_domain; db f_st5+f_two, "pow", 0

err_pow_pos:
  fld st(0), st(6)
  ftst st(0)
  fstsw [bp][-2]
  fstp st(0),st(0)
  test byte [bp][-1], 1
  jz err_plus_infinity
  jmp err_underflow

err_pow_neg:
  fld st(0), st(6)
  fistp word [bp][-2], st(0)
  mov ax, [bp][-2]
  call __pow_int
pown_exit:
  ret 0
  dw magic, pown_exit, err_domain, err_domain; db f_st5+f_two, "pow", 0

pow_dom_err:
  push si
  mov si, pow_dom_err_data
  call __lib_error
  pop si
  jmp pow_exit
pow_dom_err_data: dw err_domain, err_domain; db 0, "pow", 0


segment _TEXT(CODE,28H)
public _tanhl :
public _tanh :
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,2
tanhCommon:
  fst st(5), st(0)            (* save for error processing *)
  fldl2e st(0)
  fmulp st(1), st(0)
  fist word [bp][-2], st(0) (* to check that it is not too large *)
  fwait
  call pow2
  fld   st(0), st(4) (* retrieve integral part *)
  fxch  st(1), st(0)
  fscale st(0), st(1)
  fstp st(1), st(0)
  fld1 st(0)
  fdiv st(0), st(1)  (* Exp (-x) *)
  fst  st(7), st(0)
  fadd st(0), st(1)
  fxch st(0), st(1)
  fsub st(0), st(7)
  fdivrp st(1), st(0)
tanh_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, tanh_exit, tanh_one, tanh_minus_one; db f_st5, "tanh", 0

tanh_one:
  sub ax,ax
  fld1 st(0)
  ret 0

tanh_minus_one:
  call tanh_one
  fchs st(0)
  ret 0

segment _TEXT(CODE,28H)
public _sinhl :
public _sinh :
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,2
sinhCommon:
  fst st(5), st(0)            (* save for error processing *)
  fldl2e st(0)
  fmulp st(1), st(0)
  fist word [bp][-2], st(0) (* to check that it is not too large *)
  fwait
  call pow2
  fld   st(0), st(4) (* retrieve integral part *)
  fxch  st(1), st(0)
  fscale st(0), st(1)
  fstp st(1), st(0)
  fld1 st(0)
  fchs st(0)
  fdiv st(0), st(1)  (* - Exp (-x) *)
  faddp st(1), st(0)
  fidiv  word st(0), cs:[two]
sinh_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, sinh_exit, err_minus_infinity, err_plus_infinity; db f_st5, "sinh", 0

segment _TEXT(CODE,28H)
public _coshl :
public _cosh :
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,2
coshCommon:
  fst st(5), st(0)            (* save for error processing *)
  fabs st(0)
  fldl2e st(0)
  fmulp st(1), st(0)
  fist word [bp][-2], st(0) (* to check that it is not too large *)
  fwait
  call pow2
  fld   st(0), st(4) (* retrieve integral part *)
  fxch  st(1), st(0)
  fscale st(0), st(1)
  fstp st(1), st(0)
  fld1 st(0)
  fdiv st(0), st(1)  (* Exp (-x) *)
  faddp st(1), st(0)
  fidiv  word st(0), cs:[two]
cosh_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, cosh_exit, err_plus_infinity, err_plus_infinity; db f_st5, "cosh", 0

segment _TEXT(CODE,28H)
public _asinl :
public _asin:
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,2
asinCommon:
  fst st(5), st(0)
  fmul st(0), st(0)
  fld1 st(0)
  fsubrp st(1),st(0)
  fsqrt st(0)
  fwait
  fdivr st(0),st(5)
  fwait
  call atan0
asin_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, asin_exit, err_asin_neg, err_asin_pos; db f_st5, "asin", 0

err_asin_pos:
  call check_abs_one
  fld1 st(0)
  fstp st(1),st(0)
  fchs st(0)
  fldpi st(0)
  fscale st(0),st(1)
  sub ax,ax
  ret 0

err_asin_neg:
  call err_asin_pos
  fchs st(0)
  ret 0

segment _TEXT(CODE,28H)
public _acosl :
public _acos :
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  sub sp,2
acosCommon:
  fst st(5), st(0)
  fmul st(0), st(0)
  fld1 st(0)
  fsubrp st(1),st(0)
  fsqrt st(0)
  fdivr st(0),st(5)
  fwait
  call atan0
  fld tbyte st(0), cs:[Pi_By_2]
  fsubrp st(1), st(0)
acos_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret  0
  dw magic, acos_exit, err_acos_neg, err_acos_pos; db f_st5, "acos", 0

err_acos_pos:
  call check_abs_one
  fldz st(0)
  sub ax,ax
  ret 0

err_acos_neg:
  call check_abs_one
  fldpi st(0)
  sub ax,ax
  ret 0


segment _TEXT(CODE,28H)
public _ldexpl :
public _ldexp:
  mov  byte ss:[InMathLib], 1
  push bp
  mov bp,sp
  push ax
  fst st(5), st(0)
  fild word st(0), [bp][-2]
ldexpCommon:
  fst st(7), st(0)
  fxch st(0),st(1)
  fscale st(0), st(1)
  fstp st(1), st(0)
ldexp_exit:
  mov sp,bp
  pop bp
  mov  byte ss:[InMathLib], 0DEH
  ret 0
  dw magic, ldexp_exit, ldexp_err_neg, ldexp_err_pos; db f_st5+f_two, "ldexp", 0

segment _TEXT(CODE,28H)
ten:dw 10

public __acvt :

(* converts ascii string of digits into packed BCD then into double *)

flag    = -2
mbuf    = -6
exp     = -8
len     = -10
old_cw  = -12
stat_w  = -14
man     = -24
locals  =  10


_PNT    = 4
_MNEG   = 1
_STD    = 8
ERANGE  = 34

    mov  byte ss:[InMathLib], 1
    push  bp
    mov   bp, sp
    sub  sp, locals
    mov  [bp][flag], ax
    mov  [bp][mbuf], bx
    mov  [bp][exp], cx
    mov  [bp][len], dx
    sub   sp, 4
    mov   cx, 5
    sub   ax, ax
clearLoop:                           (* clear BCD local *)
    push  ax
    loop  clearLoop
    push  si

    mov   cx, [bp][len]
    mov   si, [bp][mbuf]
    lea   bx, [bp][man]
    shr   cx, 1
    jnc   acvtEven
    sub   ah, ah
    add   bx, cx
    inc   cx
    jmp   acvtOdd
acvtEven:
    jcxz  acvtZero
    dec   cx
    add   bx, cx
    inc   cx
acvtPack:                            (* pack number *)
    lodsb
    push  cx
    mov   cl, 4
    shl   al, cl
    pop   cx
    mov   ah, al
acvtOdd:
    lodsb
    and   al, 0FH
    or    al, ah
    mov   ss:[bx], al
    dec   bx
    loop  acvtPack
acvtZero:
    fild word st(0), cs:[ten]
    fstp st(1), st(0)
    mov ax, [bp][exp]
    test ax,ax
    jge positive_exp
    neg ax
    call __pow_int
    fbld  tbyte st(0), [bp][man]     (* load mantissa *)
    fxch  st(1), st(0)
    fdivp st(1), st(0)
    fwait
    jmp acvt_test_sign
    dw magic, acvt_err_exit, err_underflow, err_underflow; db 0, "acvt", 0
positive_exp:
    call __pow_int
    fbld  tbyte st(0), [bp][man]     (* load mantissa *)
    fxch  st(1), st(0)               (* necessary for overflow handling *)
    fmulp st(1), st(0)
acvt_test_sign:
    test  word [bp][flag], _MNEG
    je    acvt_exit
    fchs  st(0)                      (* change sign if negative *)
acvt_exit:
    pop   si
    mov   sp, bp
    pop   bp
    mov   byte ss:[InMathLib], 0DEH
    ret   0
    dw magic, acvt_err_exit, err_minus_infinity, err_plus_infinity; db 0, "acvt", 0
acvt_err_exit:
    test  word [bp][flag], _STD      (* if called from strtod *)
    jnz   acvt_exit                  (* set errno on over/underflow *)
    (* should restore errno here *)
    jmp acvt_exit

extrn __exc_addr

segment _TEXT(CODE,28H)

public __bcd_to_int :              (* si = bcd, di = int mantissa *)

  ret 0

segment _TEXT(CODE,28H)
segment _DATA(DATA,28H)

select _TEXT

__math_magic_end:
public __math_excep:
  (* called from floating point exception handler *)
  (* 8087 status word on stack *)

  push bp
  mov bp,sp
  sub sp,2
  jmp pop_8087
pop_8087_loop:
  fstp st(0), st(0)
pop_8087:
  fstsw [bp][-2]
  fwait
  and word [bp][-2], 3800H
  cmp word [bp][-2], 0800H
  jnz pop_8087_loop
  cmp byte [InMathLib], 1
  jne user_error
  mov byte [InMathLib], 0DEH

  (* search for handler *)
  push si
  push ax
  push bx
  mov si,[bp][4]
search_loop:
  inc si
  cmp word cs:[si], magic
  jne search_loop
  add si, 4
  mov ax, cs:[si][-2]    (* resume address *)
  mov [bp][4],ax     (* adjust local return address *)
  call __lib_error
  pop bx
  pop ax
  pop si
  mov sp,bp
  pop bp
  add sp, 2
  ret 0    (* discard ret to __sigFp *)

user_error:
  push ax
  mov  ax, [bp][4]
  mov  [__exc_addr], ax
  pop  ax
  mov sp, bp
  pop bp
  ret 0

end

